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Abstract 


The  behavior  of  a  communication  network  can  be  modeled  as  a  flow 
of  traffic  units  between  nodes  connected  by  links.  We  derive  a  node/link 
network  model  and  connect  it  to  a  fluid-like  model  of  traffic  flow.  The 
discrete  node/link  model  emphasizes  packet  queuing  and  the  flow  of  pack¬ 
ets  from  spatial  point  to  spatial  point.  The  model  assumes  that  packets 
reside  in  buffers  at  each  node,  and  are  classified  by  their  destination  and 
the  length  of  time  they  have  resided  in  the  buffer.  An  algorithm  was  cre¬ 
ated  for  packets  to  exit  the  buffer  at  each  node  according  to  their  age  and 
travel  to  the  next  node  along  a  predetermined  path  to  their  destination. 

This  algorithm  calculates  the  rate  at  which  packets  distribute  themselves 
to  the  next  link  in  the  route  to  their  destination,  assumes  a  source  of 
packets  originating  at  the  node,  and  subtracts  packets  whose  destination 
is  that  particular  node.  The  continuum  model  derived  from  this  discrete 
flow  model  leads  to  a  flow  continuity  equation.  The  continuity  equation 
describes  the  density  of  packets  as  a  function  of  time  and  space,  so  that 
we  are  able  to  predict  changes  in  global  flow  patterns  and  optimal  paths 
in  the  network.  Solutions  to  the  flow  equations  in  one  dimension  show 
that  if  the  sources  are  too  strong  or  the  flow  is  restricted,  the  packet  den¬ 
sity  grows  at  the  nearest  upstream  node.  When  the  source  strength  is 
reduced,  or  when  flow  is  restored,  the  buffered  packets  flow  at  capacity 
until  the  density  has  been  reduced. 

Keywords:  Communication  Networks,  Flow  Model,  Packets,  Flux,  Continuity 
Equation 


1  Introduction 

Rapid  communication  is  essential  in  many  areas  of  endeavor.  In  everyday  work 
and  pleasure,  we  communicate  important  financial  and  personal  information;  in 
military  situations,  security  and  lives  depend  on  it.  Existence  and  functioning 
in  the  information  age  requires  rapid  and  accurate  transmission  of  data  in  wire, 
optical  and  wireless  networks.  Many  systems  can  be  modeled  as  networks  where 
traffic  units  flow  along  links  connecting  nodes  [24,  21,  22].  Modeling  an  aspect  of 
flow  on  communication  networks  allows  us  to  better  understand  the  underlying 
dynamics  of  such  networks,  and  particularly  how  the  routing  of  units  (packets) 
affects  the  gross  features  of  the  flow  (delay,  congestion)  on  the  network  [22] . 

A  communication  network  is  composed  of  links  between  nodes  in  a  network 
that  carries  information  from  its  source  to  its  destination.  Everything  that  trav¬ 
els  on  communication  networks  can  be  considered  as  “quanta”  called  packets. 
Email,  web  pages,  instant  messaging  programs,  voice  communication,  etc.,  are 
all  broken  down  into  chunks  of  information  which  travel  throughout  a  commu¬ 
nication  network. 

Packets  are  generated  and  delivered  at  host  computers.  A  router  forwards 
packets  within  a  computer  network  along  an  “optimal  path”  to  another  router 
[1] .  A  packet  might  have  to  go  through  several  routers  before  it  reaches  its  des¬ 
tination.  Communication  networks  are  called  packet-switching  because  routers 
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store  and  make  forwarding  decisions  at  the  data  and  network  layers  to  direct 
packets  through  the  network. 

In  this  model  of  network  flow,  network  traffic  is  described  mathematically 
by  a  link  flow  model  for  the  buffering  of  packets  and  packet  movement  between 
discrete  sets  of  nodes  on  links  at  discrete  times.  We  shall  use  the  discrete 
description  and  an  associated  service  area  to  derive  a  continuum  description 
for  the  density  of  packets  (packets  per  square  kilometer,  for  example)  and  flow 
(packets  per  second  per  kilometer  of  flow  boundary)  which  creates  a  fluid-like 
model.  Specifically,  we  model  the  spatial  distribution  of  packet  flow.  We  seek 
out  the  variables  deemed  most  important  that  influence  traffic  flow;  allowing  us 
to  create  a  simplified  view  of  network  communication,  and  concentrate  on  packet 
distribution  as  it  depends  on  spatial  location.  For  instance,  the  only  information 
we  need  about  a  packet  at  some  location  in  the  network  is  its  destination. 

Modeling  flow  on  a  communication  network  as  a  fluid  will  allow  us  to  (*)  de¬ 
scribe  the  spatial  structure  and  evolution  of  normal  and  congested  flow  on  large 
communication  networks;  (ii)  predict  changes  in  flow  pattern  due  to  changes  in 
the  spatial  density  of  sources  and  per-link  traffic  and  (in)  determine  the  optimal 
route  a  packet  would  take  to  get  to  its  destination  (possibly  improving  upon  the 
existing  algorithms) . 

On  some  conceptual  scale,  flow  on  a  communication  network  appears  to  be¬ 
have  like  a  fluid  flow.  We  propose  a  model  fluid  flow  in  terms  of  density  of 
packets  as  a  function  of  spatial  coordinates  and  time.  We  suggest  that  this 
approach  allows  study  of  gross  features  of  flow  in  the  network  under  changes  of 
conditions.  This  report  is  organized  as  follows:  Section  2  makes  an  analogy  to 
lattice  gases  and  sets  the  structure  for  the  flow  model  that  we  wish  to  describe. 
Section  3  describes  a  buffer  model  showing  how  the  packets  evolve  in  age  and 
transition  from  buffered  to  the  send  queue.  Section  4  describes  the  discrete  com¬ 
munication  network  model.  Section  5  explains  the  continuum  model  and  how 
it  was  derived  from  the  discrete  model  which  leads  to  the  continuity  equation. 
Section  6  presents  several  solutions  to  the  model  illustrating  how  it  behaves  in 
situations  of  blockage  and  saturation. 

2  Lattice  Gas  Communication  Network  Model 

In  this  Section,  we  propose  a  model  for  the  spatial  distribution  of  packet  flow 
and  dynamics  using  the  concepts  of  the  Lattice  Gas  Cellular  Automata  model 
and  the  relation  of  lattice  gases  to  flow  equations.  Packet  motions  in  large 
communication  networks  behave  somewhat  like  a  lattice  gas.  In  a  network 
of  mobile  users,  wireless  communication  is  essential,  but  the  spatial  range  of 
wireless  is  limited.  Consequently,  messages  are  routed  from  node  to  node,  so 
that  flow  in  all  “directions”  is  possible. 

A  cellular  automata  model  is  a  system  composed  of  adjacent  cells  or  sites 
(usually  organized  as  a  regular  lattice)  which  evolves  in  discrete  time  steps. 
Each  cell  is  characterized  by  an  internal  state  whose  value  belongs  to  a  finite 
set  of  possible  packet  trafficking  and  buffering.  The  updating  of  these  states  is 
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made  at  discrete  times  according  to  a  local  rule  involving  only  a  neighborhood  of 
each  cell.  A  lattice  gas  cellular  automaton  has  specific  neighbor  interaction  rules 
and  is  a  model  of  a  gas.  Many  complex  physical  systems  such  as  traffic  flow 
can  be  modeled  as  cellular  automata.  Cellular  automata  models  can  provide 
quantitative  information  for  the  complex  systems  being  modeled. 

A  lattice  gas  model  uses  a  set  of  simple  “collision  rules”  to  model  the  much 
more  complex  interactions  of  gas  molecules.  For  the  modeling  of  a  network,  the 
role  of  routers  in  the  interaction  of  packet  dynamics  in  the  network  suggests 
different  collision  rules. 

We  derive  a  model  for  the  flow  on  a  network  accounting  for  the  role  of  routers 
or  similar  control  devices  in  determining  the  macroscopic  features  of  the  flow. 
We  will  proceed  in  analogy  with  lattice  gas  cellular  automata  by  describing 
the  interaction  of  particles  (packets)  arriving  at  nodes.  The  interactions  in 
communication  networks  differ  from  those  assumed  for  lattice  gases,  but  we 
shall  seek  common  features  which  lead  to  a  fluid  flow  description  in  terms  of 
spatial  coordinates  and  time. 

In  order  to  set  this  work  in  perspective,  we  shall  describe  the  lattice  gas  cellu¬ 
lar  automata  model,  including  the  continuum  approximation  and  the  derivation 
of  the  equation  for  the  momentum  flux  (pressure)  for  the  Euler  equation.  Then 
we  discuss  appropriate  modifications  to  make  the  flow  resemble  flow  on  a  com¬ 
munications  network.  We  derive  a  lattice  Boltzmann  approach  to  a  lattice  gas 
residing  on  an  irregular  lattice,  and  discuss  how  this  model  would  be  modified 
to  model  network  flows.  This  modeling  effort  will  lead  to  equations  for  how 
packets  move  on  a  discrete  set  of  nodes  and  links,  moving  at  discrete  times. 

2.1  Lattice  Gas-Cellular  Automata  Model 

A  lattice  gas  is  a  fluid  model  consisting  of  the  motion  of  discrete  particles  on  a 
regular  lattice  of  nodes.  A  cellular  automaton  lattice  gas  model  assumes  that  the 
gas  has  a  finite  number  of  particles  and  the  particles  move  on  links  connecting 
nodes.  At  discrete  time  steps  the  particle  moves  through  a  link  into  the  next 
node.  If  two  or  more  particles  meet  at  a  node  they  collide.  These  collisions 
are  governed  by  a  set  of  rules.  The  set  of  rules  take  into  account  the  spatial 
symmetry  of  the  lattice  and  also  conserve  mass,  momentum  and  energy. 

2.2  Lattice  Gas 

In  the  lattice  gas,  the  particle  dynamics  are  governed  by  the  equation, 

iii( x  +  Ci,  t  +  1)  =  rii(x,  t)  +  A i  [n(x,  f)]  (2.1) 

where  rij  represents  the  number  of  particles  at  lattice  point  x  at  time  t  with 
velocity  Cj  and  A;  [n(x,  t )]  are  the  collision  operators.  The  collision  operators  are 
coded  from  a  look-up  table,  and  depend  on  the  states  of  the  nodes  surrounding 
node  x.  For  the  model  of  a  gas,  the  collision  rules  are  designed  to  conserve 
mass,  momentum  and  energy.  One  representative  collision  on  a  rectangular 
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two-dimensional  lattice  is  shown  in  Figure  2.1.  In  this  model,  the  four  velocity 
directions  (i  =  1,  2, 3, 4)  are  C\  =  i,  c%  =  j,  C3  =  —i  and  C4  =  —j.  The  collisions 
are  encoded  in  the  collision  operator. 


a) 


Figure  2.1:  A  typical  collision  on  a  rectangular  lattice,  a)  before,  and  b)  after. 

We  shall  not  consider  energy  here.  Since  the  collisions  conserve  particles  and 
momentum,  we  have 
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Y  A i  [n  (x,  t)]  =  0 

2= 1 
4 

Y Ci^1  in  ( x >  *)]  =  ° 

2= 1 

Y\ciAi[n(x,t)\  =  ° 

2= 1 

Mass  Equation 

4 

Summing  the  particle  dynamics  equation  (2.1)  gives  E  Hi  ( X  +  C-ii  t  +  1)  — 


Uj  (x,£).  Multiply  the  particle  dynamic  equation  by  c*  and  sum  that  as  well 
gives 

4  4  4 

c»ni  (x  +  q,  t  +  1)  =  CiUi  (x,  £).  Define  the  density  as  p  (x,  t)  =  y^  (x,  £) 
2=1  2=1  2=1 

and  set  the  average  velocity  as 


u  (x,  t) 


4 

Y,  cini  t) 

i— 1 


p{x,t) 
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we  can  then  sum  the  particle  equation  over  all  the  nodes  inside  of  some  area  V 
(see  Figure  2.2),  and  obtain 

4  /  4 

y,  y  [n,  (a:,  f  +  1)  -  rii  (x,  t)]  =  -  I  y  y  [rii  (x  +  c*,  t  +  1)  -  rii  (x,  t,  +  1)] 
xGV  i=l  VxGV  i=l 


Figure  2.2:  A  “volume”  in  a  rectangular  lattice,  with  a  normal  vector  shown  on 
the  right. 


Note  that  the  sum  on  the  right  hand  side  is  the  net  particle  flux  out  of  the 
area  (through  the  boundary).  The  sum  on  the  right  can  be  replaced  by 

4 

E  E  v  ■  arii  ( x ,  t  +  1) 

xev  i=  1 
x+CitfiV 


where  v  is  the  unit  normal  to  the  boundary  curve.  Note  that  v  •  q  =  —1,0,1 
depending  on  whether  the  velocity  is  directed  into,  along,  or  out  of  the  boundary 

curve.  The  first  term  can  be  approximated  by  y^  - — - — ,  and  the  flux  term 
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can  be  approximated  by  v  •  p(x,t)  u(x,t).  The  spatial  sums  can  be 

x€V 

x-\-Ci(£V 

replaced  by  integrals,  an  area  integral  for  the  time  derivative,  and  a  line  integral 
for  the  flux.  We  have 


>v 


dp 

at 


dA  =  - 


v  ■  puds 


IdV 


Using  Gauss’  theorem  allows  us  to  write  this  as 


/v(s+v-'»)  "  =  0 


Assuming  smoothness  of  the  integrand,  we  have  the  equation  of  conservation  of 
mass 


dp 

dt 


+  S7  ■  pu  =  0 


Momentum  Equation 

Multiply  the  particle  equation  (2.1)  by  Ci  and  sum  gives 


4  4 

EE  Ci  [n,  (a;,  t  +  1)  -  Hi  (x,  t)]  =  -  E  E  (n  ■  ci)  CiUi  (x,  t  +  1)  (2.2) 

x€V  i—  1  x€V  i—l 

x+Citf:V 


Define  the  momentum  flux  tensor  II  =  CiCjUj  (x,  t)  and  following  the 

argument  above  gives 


i—l 


dpu 

~dt 


=  vn 


The  momentum  flux  can  be  written  as  II  =  —puu  +  a  where  the  stress  a  can 

4 

be  expressed  in  microscopic  terms  as  a  =  —  cicin.;  (x,  t).  Here  q  =  Ci  —  u  is 

i—l 

the  velocity  fluctuation.  The  momentum  equation  can  be  written  as 


p^-  +  u.Vu)=V.a 

Closure 

The  behavior  of  the  fluid  is  determined  by  specifying  the  stress  cr  in  terms  of 
the  state  variables,  which  are  expressed  in  terms  of  p  and  u.  A  simple  closure  is 
to  assume  that  the  stress  is  isotropic,  and  given  by  cr,,  =  — pSij  where  p  =  p(p,  9) 
is  the  pressure,  and  6  is  the  temperature.  This  approach  requires  the  derivation 
of  an  energy  equation,  which  is  beyond  the  purpose  of  this  discussion. 

For  lattice  gases  which  are  assumed  to  resemble  a  viscous  fluid,  closure 
is  accomplished  by  assuming  the  distribution  given  by  the  particle  dynamics 
equation  is  in  steady  state  and  homogeneous  in  space.  The  equation  becomes 
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Aj  (ar,  £)J  =  0.  For  a  viscous  flow,  the  equilibrium  distribution  is  assumed 
to  depend  on  the  shear  rate.  An  appropriate  dependence  is  assumed  for  small 
shear  rate,  and  substituted  into  the  particle  dynamics  equation.  After  assuming 
invariance,  the  coefficients  are  calculated,  and  substituted  into  the  definition  of 
stress.  This  approach  is  complicated  by  the  fact  that  the  expansion  re-introduces 
terms  expressing  momentum  flux  by  mass  flux,  i.e. ,  a  term  proportional  to  wu 

[25]. 

3  A  Probabilistic  Model  for  Flow  FIFO  Buffer 

In  this  Section,  we  descibe  a  model  for  an  individual  buffer  with  unlimited  mem¬ 
ory  operating  on  a  FIFO  (first  in,  first  out)  algorithm.  In  this  model,  packets 
are  released  based  on  their  age  a.  We  consider  the  probability  P(a,  k,  r)  of  the 
buffer  having  k  packets  at  age  a,  at  time  r.  Packets  at  age  a  can  exist  in  two 
states  depending  on  whether  a  is  the  maximum  age  of  packets  in  the  buffer. 
We  refer  to  the  two  states  as  either  “buffered”  or  “send  queue”.  The  buffered 
state  refers  to  packets  that  are  not  being  sent  but  remain  in  the  buffer.  The 
send  queue  state  refers  to  packets  at  maximum  age  in  the  buffer,  which  are  in 
the  process  of  being  sent.  Figure  3.1  shows  an  example  of  packet  behavior  as 
they  transition  from  the  buffer  to  the  send  queue.  Incoming  bunches  of  pack¬ 
ets  are  buffered  at  age  zero,  which  is  displayed  by  the  figure  on  the  right.  At 
a  certain  age  that  is  determined  probabilistically,  they  transition  to  the  send 
queue  they  begin  to  exit.  Those  packets  that  did  not  exit  age  in  the  buffer  each 
time.  The  total  probability  of  P(a,  fc,r)  of  k  packets  of  age  a  in  the  buffer  at 
time  t  is  the  probability  of  buffered  k  packets  at  age  a,  B(a ,  k.  r)  in  the  queue 
and  the  probability  of  packets  k  in  the  send  queue.  Note  that  packets  in  the 
send  queue  at  maximum  age  a  =  amax  in  the  process  of  being  sent  S(a ,  k,  r),  so 
P(a,  k,  t)  =  S(a,  k ,  r)  +  B(a ,  k,  t).  P(a ,  k,  r)  is  derived  by  understanding  each 
probability  P(a,  fc,r)  and  S(a,k,r). 

B(a ,  k,  t)  is  the  probability  of  k  packets  buffered  at  age  a.  It  is  described  by 
the  following  equation: 

B(a,  k,T  +  At)  =  B(a  —  1,  k,  r)TB(a  —  1,  k,  r)  (3-1) 

This  equation  states  that  if  k  packets  are  buffered  at  age  a  —  1  at  time  r, 
they  will  be  buffered  at  age  a,  at  time  r  +  At  with  probability  TB(a  —  1,  k,  r). 
Thus  TB(a  —  1  ,k,  r)  is  the  probability  of  transition  of  a  buffered  packet  at  age 
a  —  1  to  age  a,  implying  age  a  did  not  become  the  maximum  age  in  the  buffer. 

S(a,k,r)  is  the  probability  of  k  packets  in  the  sending  queue  at  a  =  amax 
with  the  following  equation  describing  it: 


S(a ,  k,  t+At)  =  S(a— 1,  fc+1,  r)TS'(a—  1,  k+l,  t)+B(ci—  1,  k ,  r)(l— TB{a— 1,  k,  r)) 

(3.2) 


Equation  (3.2)  states  that  the  probability  of  k  packets  going  to  the  send 
queue  is  that  either  at  time  r  they  were  already  in  the  send  queue  but  the 
number  exceeded  the  link  packet  transport  capacity,  and  therefore  were  not 
able  to  exit;  or  packets  were  previously  in  the  buffered  state  and  are  now  in 
the  send  queue.  S(a ,  k.  t)  has  two  transition  probabilities,  TS(a  —  1,  k  +  1,  r) 
and  1  —  TB(a  —  1,  k,  t).  TS(a  —  1,  k  +  1,  r)  is  the  transition  probability  that  if 
the  number  of  packets  k  is  greater  then  the  link  packet  transport  capacity,  then 
packets  that  did  not  exit  remain  in  the  buffer  and  age,  increasing  amax  by  1. 
The  probability  TS(a  —  1,  k  +  1,  r)  equals  1. 
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Figure  3.1:  Example  of  Packets  Being  in  the  Send  Queue 


The  other  transition  probability  1  —  TB(a  —  1,  k,  r)  corresponds  to  buffered 
packets  entering  the  send  queue.  This  quantity  is  the  probability  that  at  time 
r,  the  number  of  packets  at  maximum  age  amax  was  below  the  link  capacity 
(all  packets  were  sent)  and  the  maximum  age  dropped  to  a  new  lower  value 
o-max{T  +  At)  <  amaa,(r)  where  a  new  set  of  packets  that  were  previously 
buffered  will  enter  the  send  queue.  In  order  for  the  maximum  age  to  drop 
to  a  new  lower  value  there  could  not  have  been  packets  between  amaa;(r)  and 
dmax{T  +  At).  If  B(a,k,r)  is  the  probability  of  k  buffered  packets  at  age  a, 
then  1  —  B(a,k,r)  is  the  probability  of  there  being  no  buffered  packets  at  age 
a.  The  transition  probability  that  the  number  of  packets  was  below  the  link 
transport  capacity  and  the  maximum  age  dropped  to  a  new  lower  value  a  is, 


1  —  TB(a  —  1,  k,  t) 


oo  oo  a—1 

E  S(a',l,r)(l-  E  E  £(a*>fcV)) 

a'— a+1  k*  =  1  a*=a+ 1 
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The  transition  probability  TB(a  —  1  ,k,r)  in  equation  (3.1)  equals, 

oo  oo  a'  —  l 

TB(a  —  1,  fc,  t)  =  1  —  ^  5(a,,l,t)(l-53  ^  B(a*,k*,r))  (3.3) 

a/— a+1  k—1  a*  —a-\-l 

Substitute  1  —  TB(a —  l,k,r)  into  equations  (3.1)  and  (3.2),  yielding 

oo  oo  a7  — 1 

B(a,k,T+Ar)  =  B(a-l,k,T)(l-  53  5  (a',  1>T)(1“  53  53  B(a*’k*’T))) 

a'=a+ 1  /c*  —  1  a*  —a+1 

(3-4) 

oo  oo  a7  — 1 

5(a,  fc,  t+At)  =  5(a— 1,  fc+1,  r)+B(a— 1,  k,  r)  E  S(»M.0(1-E  E  B(a-,k’,r)) 

a'=a-\- 1  fc*=la*=a+l 

(3.5) 

Equations  (3.5)  and  (3.4)  fully  describes  the  two  ways  k  packets  arrive  at 
age  a.  Figure  3.2  shows  an  example  of  the  evolution  of  packets  that  are  buffered 
in  a  queue.  The  figure  on  the  left  shows  that  packets  in  this  cohort  did  not  exit 
but  rather  aged  in  the  buffer. 


Number  of  Packets  Number  of  Packets 


Figure  3.2:  Example  of  Buffered  Packets 

Equations  for  S(a,k,r)  and  B{a,k1r)  fully  exhaust  the  possible  states  for 
packets,  Buffered  or  Send  Queue.  These  two  equations  describe  packets  at  age 
a  >  1.  When  packets  arrive  at  a  node  they  enter  the  system  at  age  zero.  They 
are  either  buffered  or  enter  the  send  queue.  5(0,  k,  r  +  At)  +  B{ 0,  fc,  r  +  At) 
is  the  total  probability  of  k  packets  at  age  zero.  We  assume  that  arrivals  enter 
the  Send  Queue  if  the  system  is  empty.  The  probability  of  being  at  maximum 
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age  zero  is  the  probability  of  finding  the  system  empty, 


5(0,  k,  r  +  At) 


5(0,  k,  t  +  At)  +  5( 0,  k,  r  +  At) 


oo  /  oo  of  —  1 

=  y  s'(«,-i>t)  1  -  Z  Z  50*>fcv) 

a'=l  y  fc*  — 1  a*  =1 

(3.6) 

If  packets  do  not  enter  the  send  state  but  instead  enter  the  buffered  state  at 
age  zero  then  this  probability  is, 


5(0,  k,  t  +  At) 


5(0,  k,  t  +  At)  +  5(0,  k,  r  +  At) 


oo  a'  —  1 

=  1~Z  5(a/,l,r)  I  1-  Z  ZB(a*’fc*’r) 


k*= 1  a*  =  l 


(3.7) 


3.1  Steady  State  Analysis 
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Figure  3.3:  Age  Packet  Trajectories 

We  assume  the  probabilities  are  in  steady  state  and  the  equations  become, 
for  age  a  >  1 

oo  oo  a'  —  1 

S(a,fc)  =  5(a-l,fc)(l-  Z  S(a',l)(l-Y  Z  (3-8) 

a'=a+ 1  k*  —  l  a*—a-\- 1 

oo  oo  a  —  1 

5(a,  k)  =  S(a  —  1,  k  +  1)  +  5(a  —  1,  k)  £  s(<.',i)(i-£  £ 

a'=a+l  A:*  — 1  a*  =a+l 

(3.9) 
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and  the  equations  at  age  zero  become 


wZmtrp™  (1-f,l!s(o'’r))  (310) 

S(oZmt)  = 1  -  ?  *<“'■ » (*  ■ -  £  % <311» 

CL  —  1  \  k  —  1  CL  —  1  J 

3.2  Solution  Method 

An  approach  to  finding  the  solution  to  the  buffered  and  send  probability  equa¬ 
tions  is  by  using  separation  of  variables  for  S(a,  k)  and  B(a,  k), 

S(a,  k)  =  N(a)U(k) 

B(a,k)  =  A(a)T(k) 

Set  U(k+1)  =  XU(k)  =  X T(k).  Solving  for  C/ (At)  we  get  U{k)  =  jXk  where  we 
set  “  =  1.  This  implies  that  T(k )  =  Xk.  Substituting  into  equations  (3.8)  and 
(3.9)  yields 


oo  oo  a'  —  1 

N{a)-XN{a-l)  =  A{a-l){  -  E  E  A(a*)T(k*))) 

a'— a+1  k*—l  a*  — a+1 

(3.12) 

oo  oo  a— 1 

A(a)-A(a-1)  =  -A(a-1)(  E  N(a')U(  1)(1  -  E  E  A(a*)T(k*))) 

a/— a+1  k* =1  a*  =a-\-l 

(3.13) 

The  arrival  parameter  A  is  determined  by  the  traffic  flow  upstream  the  node 
being  described,  and  is  given  by  the  network.  Let  £7(1)  =  A  and  the  sum 
involving  T(k)  becomes  by  setting  k  =  s  +  1, 

OO  OO  OO  « 

£2W  =  5>*  =  a  £a*  =  — 

fc= 1  k—1  s=0 

Equations  (3.12)  and  (3.13)  become 


N(a)  —  AN (a  —  1) 


AA(a-l)  E 


ia/=a+l 


A(a)  —  A(a  —  1) 


— AA(a 


1) 


OO 

E 

a'=a+ 1 


A 

1  —  A 


E  ^(o 


a*=a+l 


A 

1  —  A 


E  ^(o 


a*=a+l 


(3.14) 

(3.15) 
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We  also  apply  separation  of  variables  to  equations  (3.10)  and  (3.11),  resulting 


in 


N(  0) 


N(  0)  +  A(0) 

m 

N(  0)  +  A(0) 


A 

1  —  A 
A 


a'  —  l 


a'—l 

E  ^(O) 

a*  =  1 

(3.16) 

a'  —  l 

a*  —  l 

(3.17) 

The  marginal  probability  of  age  a  is 

OO 

p{a )  =  ^(-4(a)A  +  N(a))Xk 
fc= o 

_  A(a)A  +  N(a ) 

"  1  -  A 

and  the  sum  of  all  the  marginal  probabilities  is  equal  to  1.  This  implies 

A(a)X  +  N(a) 


a—0 


a— 0 
oo 


1- A 


)+Ei4tw<“) 


a—0 

=  l 


a—0 


4  Discrete  Network  Model 

A  large  communication  network  is  composed  of  many  hundreds  to  hundreds  of 
thousands  of  routers  distributing  packets  throughout  the  system.  Graph  theory 
describes  the  connectivity  of  the  network  [23,  24].  A  graph  is  composed  of  a  set 
N  of  nodes  (vertices)  and  a  set  E  of  edges  (links)  where  an  edge  connects  a  pair 
of  nodes.  In  consequence,  G  =  (N,  E)  is  the  underlying  structure  of  a  network. 
In  a  communication  network,  information  travels  between  the  nodes  along  the 
links  that  connect  them.  These  links  may  be  directed  or  undirected.  In  directed 
graphs  a  link  has  a  preferred  direction  for  information  to  travel,  while  in  an 
undirected  graph  there  is  no  distinction.  The  graph  in  this  model  is  undirected; 
information  travels  in  both  directions.  In  this  network  model  a  matrix,  called 
a  route  matrix,  describes  the  path  from  source  to  destination  for  a  packet  to 
travel  optimally  in  the  network. 

We  simplify  the  model  by  assuming  that  routers  are  both  host  and  switch 
computers,  implying  that  packets  are  created,  destroyed  and  buffered  at  the 
routers.  Packets  will  flow  through  links  to  each  router.  Upon  entering  the 
network,  a  packet  will  go  directly  into  their  source  buffer  and  will  eventually 
exit  the  network  at  a  pre-determined  router.  Each  router  instructs  each  resident 
packet  which  link  to  take  as  it  travels  to  its  destination. 
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4.1  Buffer  State 


The  flow  on  this  network  is  described  by  specifying  b{j1  d,  a,  r)  the  number  of 
buffered  packets  at  any  time  r  at  each  node  j  £  N,  with  destination  d  £  N,  of 
age  a  time  units  spent  in  the  buffer. 

4.2  Route  Matrix 

The  network  studied  here  is  assumed  to  have  packet  routes  determined  statically. 
The  route  that  a  packet  takes  from  router  j  to  its  destination  d  is  determined 
based  on  predetermined  knowledge  of  the  network.  We  do  not  treat  algorithms 
that  change  routes  dynamically  in  response  to  traffic  patterns.  Specifically  when 
a  packet  enters  a  network  at  a  router,  its  complete  path  to  its  destination  is 
determined. 

The  route  matrix  Rjtd  is  an  N  x  N  matrix,  where  N  is  the  total  number  of 
routers  in  the  network.  The  entries  Rj,d  of  the  route  matrix  are  the  appropriate 
next  router  a  packet  will  take  as  it  propagates  toward  its  destination.  Specifi¬ 
cally,  Rhd  =  x  implies  that  for  a  packet  at  router  j  with  address  (destination) 
d,  x  is  the  next  router  in  the  path  [27]. 

4.3  Queue  Dynamics 

Each  router  contains  a  queue  of  buffered  packets.  Packets  are  generated  and 
delivered  at  each  router.  When  a  packet  is  generated  it  is  given  a  destination 
by  its  originator.  Packets  are  assumed  to  have  the  same  size.  Packets  enter  a 
queue  in  one  of  two  ways: 

1.  Flowing  from  another  router  along  a  connecting  link. 

2.  Being  generated  at  the  router  according  to  some  given  packet  generator 
rate  appropriate  to  that  router. 

Once  a  packet  reaches  its  destination  it  is  delivered,  and  exits  immediately. 
We  expect  the  number  of  packets  in  a  queue  to  remain  relatively  low;  otherwise, 
it  is  an  indication  of  congestion.  The  number  of  packets  entering  the  buffer  at 
router  j  is  described  by  the  following  equation, 

oo  N 

b(j,  d,  0,  r  +  1)  =  ^2  ^2  d,  a, T )  +  d ,  r)  (4.1) 

a— 0  i=i 

Here  (3{i ,  d,  a ,  t)  denotes  the  number  of  packets  sent  from  the  node  i  toward 
destination  d  at  time  r.  Each  packet  sent  is  chosen  from  the  population  of 
packets  at  the  upstream  nodes  characterized  by  having  age  a  in  that  buffer.  The 
Kronecker  delta  connects  the  next  location  j  from  node  i  with  destination  d  in 
terms  of  the  route  matrix.  The  second  term  on  the  right  side,  u(j,  d ,  r),  describes 
the  number  of  new  packets  entering  the  system  at  node  j,  with  destination  d, 
at  time  r.  Packets  entering  the  buffer  have  age  a  =  0. 
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We  assume  that  each  router  uses  a  FIFO  unlimited  memory  buffer  where 
packets  are  released  based  on  their  age  a.  The  number  of  packets  in  the  buffer 
are  described  by  the  evolution  equation 

K h  d,  a  +  1,  r  +  1)  =  b(j,  d,  a,  r)  -  0(j,  d,  a,  r)  (4.2) 

The  second  term  on  the  right  hand  side  j3 (j,  d,  a ,  r),  determines  the  number 
of  packets  that  are  sent  outward  from  node  j  with  destination  d  at  each  time. 

At  any  router,  destinations  are  competing  for  the  same  link  to  a  router.  This 
creates  competition  among  destinations  as  to  the  number  of  packets  that  will 
have  access  to  the  link.  Due  to  this  competition  an  algorithm  was  created  to 
describe  the  flow  of  packets  to  a  link.  We  note  that  if  a  packet  buffered  at  j 
with  destination  d  shares  the  next  link  Rj  d  with  a  packet  buffered  at  j  with 
destination  l,  there  is  a  possible  conflict  with  sending  both  packets. 

4.4  The  Flow  Variable 

To  determine  the  number  of  packets  that  flow  /3(j,  d,  a,  r)  from  node  j,  with 
destination  d,  we  have  to  determine  the  number  of  packets  that  are  competing 
at  each  time  for  a  link.  We  define  a  function  k(j,  d,  a ,  r),  which  sums  the  total 
number  of  packets  having  destinations  which  compete  for  the  same  next  link  at 
age  a.  The  total  number  of  packets  k(j.  d,  a,  r),  to  be  sent  from  j  to  Rh,i  at  age 
a  is, 


N 

k(j,  d,  a,r)=Y^  Kj,  h  a,  T)6Rj  ltRj  d  (4.3) 

l=i 

The  Kronecker  delta  in  equation  (4.3)  separates  the  links  that  are  competing 
with  destination  d  by  using  the  route  matrix.  The  buffer  at  a  router  needs  to 
know  which  destinations  are  competing  for  the  same  link  to  the  next  router, 
since  each  link  has  a  packet  capacity.  If  destinations  are  competing  and  there 
are  packets  ready  to  exit,  we  use  Kronecker  delta  function  to  distinguish  which 
destinations  are  competing  for  the  same  next  router. 

Mathematically,  the  competition  for  two  destinations  d  and  l  is  given  by 


1,  if  Rj,d  .1  ■ 

0,  if  Rj,d  ^  Rj,i  ■ 


The  number  of  packets  flowing  from  node  j  to  node  Rj.d  is  modeled  as, 


0{j,  d ,  a,  t)  =  b(j,  d,  a,  r)  c(j,  d,  a,  r)  (4.4) 

j3{j,  d,  cl,t)  is  the  number  of  packets  at  age  a  with  destination  d  that  will  be 
flowing  away  from  node  j.  The  competition  function  c(j,  d ,  a ,  r)  determines  the 
fraction  of  packets  at  age  a  that  will  be  sent.  The  competition  function  equals 


c(j,  d,  a,  t ) 


bo, 

bo  -  K(j,  d,a,r) 
k(j,  d,  a,  t) 


0, 
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if  a  >  a*; 
if  a  =  a*; 
if  a  <  a*. 


where  bo  is  the  bandwidth  or  total  link  packet  transport  capacity.  The  function 
k(j ,  d ,  a,  t)  determines  the  total  number  of  packets  at  each  age  a  competing  for 
the  same  next  link,  a*  is  the  minimum  age  in  the  buffer  that  can  send  packets, 
and 

OO 

K(j,  d,a,r)  =  k(j,  d,d,r) 

a=a-\- 1 

is  the  number  of  packets  above  age  a.  Consequently,  if  for  some  a,  K (j,  d ,  a,  r)  < 
bo  and  K(j,  d,  a  —  1,  r)  >  bo,  then  we  define  a*  =  a.  All  packets  at  ages  greater 
than  a*  are  sent,  hence  c(j,d,a,r)  =  1  for  a  >  a* .  However,  if  the  number  of 
packets  being  sent  at  ages  above  a*  is  below  b0  then  some  of  the  packets  at  a* 
also  contribute  to  the  total  number  of  packets  that  can  be  sent  up  to  the  capacity 
bo-  Then  b0  —  K(j ,  d,  a* ,  r)  is  the  excess  number  or  packets  that  can  be  sent  at 
age  a*,  and  since  destinations  are  competing  for  the  same  next  link,  c(j ,  d,  a* ,  r) 
determines  the  fraction  of  packets  that  can  be  sent  at  each  destination.  Packets 
with  age  less  than  a*  will  not  be  sent  since  the  total  number  of  packets  being 
sent  has  reached  link  capacity  and  consequently  /3(j,d,a,r)  =  0  for  a  <  a*.  If 
a*  =  0,  then  K(j,  d,  0,  r)  <  b0  resulting  in  all  packets  exiting. 

4.5  Discrete  Flow  Equation 

The  number  of  packets  at  router  j  evolves  in  time  as  packets  arc  flowing  through¬ 
out  a  network.  Let  n(j,  d ,  r)  be  the  number  of  buffered  messages  at  router  j 
with  destination  d  at  time  r.  Then 


OO 

n(j,  d,r)  =  a’  r)  (4.5) 

a— 0 

so  that  the  number  of  packets  at  router  j  changes  in  time  according  to 

N  N 

n{j,  d,r+  1)  —  n{j ,  d,  r)  =  ^  b{j,  d,  a,  r  +  1)  -  ^  b(j,  d ,  a,  r)  (4.6) 

a— 0  a— 0 

Using  the  evolution  equations  (4.2),  equation  (4.6)  can  be  rewritten  and  ex¬ 
panded  as 


n(j,  d,  t  +  1)  -  n(j,  d,  t)  =  ^(&(j,d,a,r)  -  f3(j,d,a,r)) 

a— 0 

oo 

+b(j,  d,  0,  t  +  1)  —  ^2  Kj,  d,  a,  r)  (4.7) 

a— 0 

Also  substituting  equation  (4.1)  for  b(j,  d,  0,r  +  1)  into  equation  (4.7),  we 
have 
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oo  oo  N 

n(j,d,T  +  1)  —  n(j,d,T)  =  -  ^  /3(j,  d,  a,  t)  +  ^  ^  SjiRi  df3{i,  d,  a,  r) 

a— 0  a— 0  i=i 

OO 

+^(i,  d,T)  •?’ a’  r)  (4-8) 

a— 0 

Equation  (4.8)  is  the  discrete  equation  expressing  “conservation”  of  packets 
at  node  j.  The  first  term  on  the  right  side  of  the  equal  sign  is  the  number  of 
packets  that  flowed  away  from  node  j.  The  second  term  describes  the  number 
of  new  packets  into  node  j  from  the  surrounding  connected  nodes.  Packets 
whose  destination  is  router  j  exit  immediately.  This  is  expressed  by  the  term, 
—(d{j-  j,  a.  t )  which  becomes  the  sink  term  at  time  r.  The  third  term  is  the 
source. 


5  Continuum  Description 

5.1  Geometry 

The  spatial  location  of  each  router  in  a  computer  network  is  of  importance 
in  some  applications.  A  router  serves  a  particular  coverage  area,  distributing 
packets  to  and  from  users  who  are  using  the  network.  To  describe  the  spatial 
and  temporal  behavior  of  a  network,  routers  are  associated  with  physical  spatial 
locations.  Let  xj  be  the  location  of  router  j.  Let  y  be  the  location  of  destination 
d.  Each  Xj  generates  a  Voronoi  polygon  V(x.,'),  a  tessellation  of  the  plane  [16]. 
The  Voronoi  polygon  V(xj)  is  the  set  of  points  closer  to  the  point  (router)  Xj, 
than  to  all  the  other  routers  with  location  x/.,  k  ^  j  in  the  domain. 

Each  router  Xj  shares  physical  boundaries  with  another  set  of  routers  x&(xj) 
so  that  packets  moving  through  the  network  will  pass  through  physical  bound¬ 
aries  around  routers.  Points  on  the  boundary  of  the  Voronoi  polygon  are  equal 
distance  from  the  two  routers  that  define  the  boundary. 

Given  P,  a  finite  set  of  distinct  points  in  the  Euclidean  plane,  each  point 
generates  a  Voronoi  polygon  in  the  plane.  A  collection  of  such  polygons  in 
Euclidean  space  define  a  Voronoi  Diagram  Vp.  The  set  given  by 

Vp  =  U  V(*j) 

j£P 

is  the  planar  ordinary  Voronoi  diagram  generated  by  P  (or  the  Voronoi  diagram 
of  P).  We  denote  the  area  of  the  Voronoi  polygon  by  |  V (xj)  |. 

When  a  network  is  created  a  route  matrix  is  determined  for  optimal  paths. 
The  route  matrix  indicates  which  routers  are  communicating  directly  with  each 
other.  Communicating  routers  do  not  necessarily  share  a  boundary.  Bound¬ 
aries  are  defined  by  the  Voronoi  polygons  generated  by  a  router  and  its  spatial 
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neighbors.  We  are  interested  in  the  rate  of  change  of  the  number  of  pack¬ 
ets  in  a  Voronoi  diagram  Vp  and  the  flow  of  packets  entering  and  exiting  the 
boundaries.  We  are  not,  however,  interested  in  the  flow  internal  to  the  Voronoi 
diagram  between  routers  and  their  Voronoi  polygons  V(xj),j  G  P. 

In  order  to  describe  the  flow  into  and  out  of  a  spatial  region,  consider  a 
Voronoi  diagram  Vp  and  its  boundary  dVp.  The  boundary  is  made  up  of  several 
straight  line  segments,  each  of  which  is  a  shared  boundary  between  the  two 
Voronoi  polygons,  V(xi)  and  V(xm)  with  l  G  P  and  to  ^  P.  We  denote  the 
boundary  element  by  dVi/rn,  and  its  length  by  |  dVi>m  |.  Then  dVp  =  [J  dVi.m 

l£P 

m(£P 


5.2  Density 

Routers  are  endowed  with  physical  spatial  locations  in  order  to  derive  flow 
equations.  Let  Vp  be  a  Voronoi  diagram  and  d Vp  its  boundary.  If  n(j,  d,  r) 
is  the  number  of  buffered  packets  at  node  j  (location  Xj)  with  destination  d 
(location  y),  then  we  define  p(xj,y,t)  as  the  density  of  packets  associated  with 
Voronoi  polygon  V  (x7)  by 


P(xj,  y,i) 


n(j,  d,  t ) 
I  V(xj)  | 


(5.1) 


where  |  V (x7j  |  denotes  the  area  of  V  (xj)  and  t  =  rAt  and  0  <  t  <  oo. 

We  assume  that  p(xj,y,t)  is  a  smooth  function  of  the  spatial  variables  Xj 
and  y  and  time  t.  Consider  the  number  of  packets  buffered  in  the  Voronoi 
diagram  Vp: 

^2  n(j ,  d,  t) 

jeP 


which  then  equals 


nC?>  d’T)  =  Y  I  V  (xj)  I 

jep  j&p 

The  continuum  limit  of  this  quantity  is 


^  n(j,  d, r) /  p(x,y,t)dV  (5.2) 

jeP  Jvp 

where  — »  denotes  an  appropriate  continuum  limit. 


5.3  Flow  Equation 

We  derive  an  evolution  equation  for  p(x;,  y,  t)  by  summing  equation  (4.8)  over 
j  G  P  and  dividing  by  At.  Thus, 


'  Vp 


P(x,  y,  t  +  At)-  p(x,  y,  t)  dv  =  ^2  n(p  d,T+l)  -  n(j,  d,  t)  ^  ^ 


At 


jeP 
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where 


E 

j£P 


n(j,  d,  t  +  1)  -  n(j,  d,  r) 
At 


uu  ^ 

jGPa=0  ^ 

N 

'52S3,Ri,d0(i>d,  a>r) 

i£P 

+v(j,d,r)  -  0{j,j,a,T))  (5.4) 


If  we  take  the  limit  as  At 
becomes 


goes  to  zero,  the  first  term  in  equation  (5.4) 


5.4  Flux 


We  consider  the  flux  terms  in  the  evolution  equation.  We  focus  on  one  Voronoi 
diagram  boundary  element  dV  =  where  l  £  P  and  m  £  P  with  l  and 

to  share  a  boundary  of  their  Voronoi  polygons.  We  define  the  flow  direction 
between  routers  in  terms  of  the  unit  vectors  connecting  router  xt  to  router  x; 

X,;  —  X; 

as  , — - The  flow  vector  from  node  x,  to  node  x^  is  then 

I  x»  -  Xj  | 


7  0{j,d,a,T) 


Note  that  i  =  Rjtd,  and  define 

_  J  1,  if  the  line  segment  connecting  x;  to  x,  passes  through  9V)jm; 
*)Kz,m)  —  1  0,  otherwise. 


This  counter  determines  whether  the  flux  from  router  Xj  to  router  x,  passes 
through  Define  rp>m  as  the  unit  normal  vector  on  boundary  point¬ 

ing  in  the  direction  of  the  adjacent  router  xm  outside  of  V  but  communicating 

xm  —  X; 

with  x; .  Then  rp  m  =  - - -  is  the  unit  exterior  normal  to  boundary 

I  xm  -  X/  | 

element. 

In  equation  (5.4)  the  outgoing  flux  term  becomes 


EE 

a= 0  j£P 
i<£P 


-1  X,;  -  Xj 

At  |  X,;  -  Xj  | 


0{j,  d,  a ,  T ) 


We  can  define  the  outflow  per  unit  area  $  through  boundary  element  dVi,m 
with  length  |  dVitm  |  as 


o°  i 

'y  v  y  v  (i,m)di,Rj  d  XT  ]  r d,  a,  t) -ip)TO  =  $  (x/jm,  y,  i) ■  n ;,m  |  dVi^m  \ 

a—0  jep  ' 

it£P 

(5.5) 
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where  x/  m  is  the  midpoint  of  dVi  m,  the  unit  vector  normal  n;  m  = 

|xi  Xj  | 

and  is  outgoing  flux.  Again,  we  assume  that  4>°(x;im,  y,  f)  is  a 

smooth  function  of  x,  y  and  t. 

A  similar  equation  for  the  incoming  flux  term  is  given  by 


OO 

a— 0 

jep 


1  Xj  -  Xj 

At  |  x.j  -  xj  | 


d,  a,  r)-ni  m  = 


<1 >J(x,  ,7715  y,i)‘  ,771  I  av;  ,771  | 

(5.6) 


The  flux  of  packets  entering  and  exiting  the  boundary  <9V,m  is 
3>(x; 

,771  5  y  5  £)-n/ 

,771  I  c>V  ,771  |  —  *°(xj  ,7715  ,771  I  dVi  ,771  |  ^J(x;  ,7715  y  1  £)•  n;  ,771  I  c>V  ,771  | 

with  <f>/(xjjm,  y,  t)  as  the  incoming  flux  terms.  When  we  sum  over  all  2  €  P  and 
to  (jL  P,  we  get  the  total  flux  of  packets  entering  and  exiting  a  Voronoi  diagram 
Vp  with  destination  y 


,7715  y  1  t)  ■  n I  5V,m  I  (5.7) 

iedP 

m€dPi 

where  <9P  =  {2  6  P  V(x;)  fl  V(xm)  g  0  |  to  ^  P}  and  <9P;  =  {to  ^  P  | 

v(x;)  n  v(xm)  ^  0}. 

This  sum  gives  the  net  flux  across  the  boundary  of  a  Voronoi  Diagram  Vp 
into  Vp,  which  becomes  an  integral  in  the  continuum  limit 


iedP 

m£dPi 


ni 


dVx, 


$(x,y ,t)  -nds 


where  |  <9V,m  |  becomes  ds. 


5.5  Source  and  Sink 


The  derivation  for  the  source  and  sink  terms  are  very  similar.  We  can  sum  the 
v(j,  d ,  r)  terms,  the  source  term  for  new  packets  entering  the  Voronoi  Diagram 
Vp,  with  destination  y,  from  equation  (4.8) 


je  p 


At  |  V(xj) 


We  can  do  the  same  thing  for  the  sink  terms 

At  |  V(x,)  |  1  V[Xj) 

CL= 0  jeP  1  v  3J  1 

If  we  define  the  continuum  limit  of 

1  P(j,j,a,T) 


At  |  V(xj) 


=  7(x,y  ,t) 


(5.8) 


(5.9) 


(5.10) 
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we  have 


(5.11) 


E  &tv{j'd'T) 


j£P 


/  Vp 


7(x,y,t)  dv 


We  can  do  the  same  to  the  sink  term,  resulting  in 


EE 


a— 0  j£P 


/  Vp 


r(x,  t)  (A 


(5.12) 


5.6  Continuity  Equation 

Combining  the  integral  of  the  density  changing  in  time,  the  flux  term  (5.7) 
and  the  source  and  sink  terms  (5.11)  and  (5.12)  we  arrive  at  the  equation  for 
conservation  of  packets  for  Voronoi  Diagram  Vp: 


4/  P(x,  y,  t)  dV  =  -  [  $(x,y,f)-nds+  f  j(x,  y,  t)  dV  -  f  a(x,t)dV 

ut  JVp  J dVp  JVp  JVp 

(5-13) 

Assuming  sufficient  smoothness  and  using  Leibnitz  rule  and  the  divergence  the¬ 
orem  gives 


J  +  v  •  $(x,y,f)  -7(x,y,t)  +  a(x,t))dV  =  0  (5.14) 

Again,  for  sufficient  smoothness  of  the  density,  flux,  source  density  and  sink 
density,  we  can  use  the  Dubois  Reymond-  Lemma  to  derive  the  continuity  equa¬ 
tion, 


+  V  •  3>(x,  y,  t)  =  q(x,  y,  t)  -  a(x,  t)  (5.15) 

5.7  Connection:  Transport  Model  and  Cellular  Automata 
Lattice  Gas  Model 

When  we  consider  a  network  of  nodes  connected  by  links,  the  spatial  location 
of  the  nodes  is  not  regular,  and  the  links  are  not  restricted  to  lie  in  regular 
directions.  In  analogy  with  the  cellular  automata  model  for  a  gas,  we  shall 
assume  that  packets  are  of  uniform  size,  and  are  not  buffered  but  sent  out 
immediately.  In  this  case,  b(j,  d,  0,  r)  =  n(j,  d,  r)  and  b(j,  d ,  a,  r)  =  0  for  a  >  0. 
Moreover,  the  velocity  of  the  particles  (packets)  is  determined  by  the  physical 
locations  and  the  propagation  speed  of  the  links. 

In  analogy  with  the  particle  dynamic  equation  (2.1),  we  apply  the  above 
assumptions  to  equation  (4.8)  without  sources  or  sinks.  We  also  note  that  in 
this  case,  the  age  variable  a  is  always  equal  to  zero,  since  no  packets  are  buffered, 
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and  the  flow  variable  is 


oo  oo  N 

n(j,  d,r  +  l)  —  n(j,  d,  r )  =  -^2 /3(j,d,a,r)  +  ’Y^2,5jtRi  df3(i,d,a,T) 

a— 0  a— 0  i= l 

j^d 

N 

=  ~bo b(j,  d ,  0,  t)  +  ^2  t>j,Ri,dbob(i,  d,  0,  r) 

J^d 

N 

=  ~bon(j,  d1r)  +  '22  Sj,Ri,dbon{i,  d,  r) 

Figure  (5.1)  shows  the  “collision”  model.  Continuing  the  analogy  with  the 
lattice  gas,  we  define  the  velocity  along  each  link  as  c,.,  =  bo  |x‘-xJ|  >  an<^  eflua_ 
tion  (5.16)  can  be  written  as 

N 

n(j,d,r+  1)  -  n(j,  d ,  r)  -  Cj,R.  dn(j,  d,  r)  •  ^  +  ^r.  dCi jn(i,  d,  r) 

j?d 


Figure  5.1:  Typical  interaction  on  a  non-uniform  lattice. 


6  One-Dimensional  Flow  Model 

Let  us  consider  the  one-dimensional  continuity  equation  flow  between  a  given 
source,  taken  as  x  =  0  and  with  a  single  destination  at  x  =  y.  The  equation  is 
then 

%  +  §xHx,v,t)  =  o  (6.i) 

for  0  <  x  <  y,  with  p(x,  y,t)  =  0  otherwise. 
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If  we  assume  that  y  >  x,  the  flow  goes  from  left  to  right  and  the  left  and  right 
boundary  points  are  the  source  and  sink  respectively.  We  are  interested  how  the 
density  changes  at  the  nodes  between  the  source  and  the  sink.  The  equations 
are  non-dimensionalized  so  that  in  the  non-saturated  one  dimensional  model 
the  flow  variable  j3(j,a,r)  is  equal  to  the  density  p(x,t)  and  the  maximum 
saturation  flow  is  equal  to  one. 

We  shall  assume  that  the  flux  in  one  dimension  in  this  model  is  given  by 

f\  =  f  p(x,  y,  t)  if  p(x,  y,  t )  <  <f>max(x,  y ,  t) 

y  \  ®max{x,y,t)  if  p(x,y,t)  >  $max(x,y,t) 

The  flow  saturation  is  $mra(a;,  y,t)  is  assumed  to  be  <  1.  We  refer  to 
&max(x,y,t)  =  1  as  system  saturation.  We  wish  to  study  situations  where 
®max(x,y,t)  <  1  for  limited  regions  and/or  time  intervals,  representing 

We  will  use  the  method  of  characteristics  to  solve  the  continuity  equation. 
When  the  flow  is  non-saturated,  then  $(x,  t)  =  p(x,t).  Our  one-dimensional 
continuity  equation  (20)  becomes 


dp 

dt 


(6.2) 


When  we  solve  equation  (21)  using  the  method  of  characteristics  and  with  initial 
condition  p$(x,Q),  our  solution  to  equation  (21)  is  p(x,t)  =  po{x  —  t).  This 
solution  is  valid  unless  there  is  a  discontinuity  in  the  density  or  the  flux.  We 
explore  several  cases  below. 


6.1  Non  Saturated  Flow 

In  this  example  we  begin  with  an  initial  density  of  packets,  and  as  time  pro¬ 
gresses  no  new  packets  are  introduced  into  the  flow.  We  will  examine  how  the 
flow  is  progressing  and  watch  how  the  density  at  each  spatial  grid  point  changes. 
Since  the  flux  density  $(x,  t)  =  p(x,t),  no  packets  should  stay  at  a  spatial  grid 
point  since  they  are  all  below  the  link  packet  transport  capacity  $mal(X;  t).  The 
function  p0  is  given  by 


Po(x) 


0.2  40  <  x  <  80 

0.4  0  <  x  <  40 

0  x  <  0 


Figure  6.1  displays  the  initial  condition  and  Figures  6.2,  6.3  and  6.4  show 
the  flow  progressing  in  time  as  a  wave  moving  from  left  to  right. 

For  this  simple  flow,  the  model  performed  as  expected. 
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Figure  6.1:  Packet  Density  ,  Initial  Condition 


Figure  6.2:  Packet  Density,  Time:  0,10,20 


6.2  Flow  with  Saturation 

The  flow  model  in  the  previous  section  dealt  with  a  flow  where  p{ x,  y ,  t)  was 
less  then  the  link  packet  transport  capacity.  Using  the  same  initial  condition  as 
the  previous  section  we  will  examine  an  example  where  p(x,  y:  t)  can  be  greater 
than  the  link  packet  transport  capacity.  We  assume  that  the  source  point  x  =  0 
is  emitting  packets  into  the  system  as  follows 


v 


(0,  80,  t)  = 


1.25  t  <  30 
0.5  t  >  30 


We  now  have  the  solution 


p{x,  t.)  =  S(x)f(t)  +  p0(x  -  t) 


24 


Density 


Figure  6.3:  Packet  Density  Scenario,  Time:  20,30,40,  Time:  40,60,80 
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Figure  6.4:  Packet  Density,  Time:  80,90,100 
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Figure  6.5:  Packet  Density  Scenario,  Time:  0,10,20,  Time:  20,30,40 


Figure  6.6:  Packet  Density  Scenario,  Time:  40,60,80,  Time:  80,90,100 


Examining  the  flow  with  saturation  we  notice  that  these  are  packets  buffered 
at  x=0,  and  that  the  flow  remains  at  the  link  packet  transport  capacity  for  some 
time  for  an  interval  to  the  right  of  a:  =  0,  but  the  buffered  packets  clear  the 
node  at  x  =  0  after,  and  the  flow  eventually  remains  at  0.5. 

6.3  Interruptions:  Disturbances  in  the  Flow 

One  of  the  issues  we  wanted  to  understand  in  this  network  model  is  what  hap¬ 
pens  to  the  network  flow  when  there  is  a  disruption  of  some  sort  in  the  flow. 
In  this  example  we  examine  interruptions  to  the  flow  in  one  dimensional  mod¬ 
els.  Interruptions  in  this  model  might  occur  in  one  of  two  types,  (i)  limited 
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bandwidth  or  more  specifically,  the  link  packet  transport  capacity  drops  to  a 
lower  value;  or  (ii)  router  interruptions,  when  the  routers  corresponding  to  a 
grid  point  are  down  for  some  time.  We  use  the  same  initial  condition  as  before 
to  investigate  interruptions.  This  model  scheme  may  has  several  interruptions 
that  can  take  place. 

The  source  node  has  a  constant  source  resulting  in  a  density  downstream 
corresponding  to  0.4  packets.  Down  town  is  experienced  around  the  point  x  = 
15.  This  down  time  equals  4>maa,(;r,  t)  equals 


given  by 
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The  figures  are  side  by  side  in  this  example  where  one  shows  the  density 
in  “normal”  resolution,  while  the  other  shows  a  close  up  of  the  neighborhood 
of  the  grid  points  where  <5>max(%,  t)  =  1. 


Figure  6.7:  Packet  Density  Scenario,  Time:  0,10,20 

Figures  6.7,  6.8,  6.9,  and  6.10  clearly  demonstrate  what  is  happening  in  this 
example.  First  we  see  the  density  is  growing  at  the  inner  grid  point  x  =  15 
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Figure  6.8:  Packet  Density  Scenario,  Time:  20,30,40 


Figure  6.9:  Packet  Density  Scenario,  Time:  40,60,80 

while  it  is  down  for  a  period  of  time,  which  is  clearly  demonstrated  in  the 
figures.  When  the  grid  point  x  =  15  is  returns  to  fully  functioning  the  flow 
recovers  quickly. 
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Figure  6.10:  Packet  Density  Scenario,  Time:  80,90,100 
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7  Conclusion 


The  goal  was  to  derive  a  model  for  the  fluid  like  aspects  of  flow  in  a  communi¬ 
cation  network,  and  to  study  the  predictions  of  such  a  model  under  conditions 
of  high  flow  and  interruptions.  The  model  derived  in  this  research  modeled  the 
behavior  of  nodes  representing  routers  forwarding  packets  along  links  to  other 
routers.  A  model  for  the  behavior  of  a  single  buffer  operating  a  first  in,  first 
out  queuing  protocol  using  a  route  matrix.  Note  that  this  matrix  appears  in 
the  flux  terms.  This  requires  that  each  packet  in  the  network  is  also  identified 
by  its  destination.  As  packets  flowed  through  a  network  from  a  source  to  their 
destination,  along  the  way  they  are  buffered  at  next  appropriate  routers  along 
their  paths.  Counting  the  number  of  packets  that  enter  and  exit  the  network 
at  a  router,  as  well  as  the  number  of  packets  that  flowed  to  and  from  led  to 
the  discrete  “conservation”  of  packets  equations.  In  order  to  view  this  model 
as  a  fluid  flow  model,  spatial  location  of  the  routers  were  defined.  Defining  the 
spatial  location  of  routers  led  to  a  continuum  description  of  the  other  variables 
that  describe  network  flow  in  this  model,  namely  density,  flux  etc. 

The  communication  network  model  derived  in  this  report  describes  the  flow 
of  packets  with  a  continuity  equation,  giving  the  density  of  packets  as  a  function 
of  time  and  space.  The  concept  of  momentum  balance,  while  not  finished  in 
this  work,  can  be  applied  to  this  model.  The  analogy  with  <9ti>  +  v  ■  Vw  =  -Vw 
[25]  is  clear;  however,  it  is  less  clear  what  role  stress  and/or  pressure  have  for 
communication  network  flows. 

We  also  explored  aspects  of  the  analogy  of  network  traffic  flow  with  cellular 
automata  models.  Cellular  Automata  models  attempt  to  simplify  large  physical 
systems  by  allowing  space  and  time  to  be  discrete  and  the  physical  quantities 
to  take  only  a  finite  set  of  values.  This  approach  describes  gross  flow  patterns 
resulting  from  packet  dynamics. 

This  communication  network  model  consists  of  a  lattice  comprised  of  nodes 
and  links.  The  particles  on  the  lattice  represent  the  packets  in  the  internet. 
Interaction  of  the  packets  takes  place  through  the  logic  implemented  by  the 
routers.  Each  link  is  a  connection  between  two  nodes  and  packets  will  travel 
along  predetermined  links  to  get  to  their  destinations.  The  rules  governing  the 
interaction  of  packets  at  nodes  on  the  lattice  in  this  model  are  those  implied  by 
equation  (5.16),  so  that  “collisions”  will  be  represented  by  changes  in  directions 
implied  by  the  route  matric. 

This  work  studied  one-dimensional  numerical  models  applied  to  various  sce¬ 
narios  to  examine  the  behavior  of  packet  flows.  Although  some  interesting 
behaviors  have  been  found  with  the  one-dimensional  model  there  are  still  many 
things  that  can  be  studied.  For  instance,  instead  of  the  buffers  having  infinite 
capacity,  the  model  could  assume  that  they  are  finite  and  when  the  densities 
reaches  some  maximum  then  packets  are  lost.  Furthermore,  a  realistic  model 
might  assume  that  when  a  buffer  is  full  then  routers  from  upstream  cannot  send 
packets.  This  then  would  represent  a  more  believable  backup  in  the  system. 

All  models  simulated  the  flow  for  an  interval  0  <  x  <  80  with  0  <  t  <  100. 
In  that  domain  the  models  demonstrated  interesting  and  physically  relevant 
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behavior.  In  the  first  example  we  focused  on  flows  without  saturation,  using 
these  calculations  to  make  sure  the  model  performed  properly. 

The  second  example  demonstrated  the  flow  with  saturation.  These  simula¬ 
tions  demonstrate  that  the  maximum  flux  model  seems  to  behave  in  a  reasonable 
manner.  In  the  last  example  we  studied  the  effects  of  interruptions.  Indeed  we 
created  network  flow  model  to  study  the  effects  of  interruptions  in  a  network.  In 
these  one-dimensional  models,  interruptions  had  a  predictable  effect  on  the  flow, 
causing  backups  in  the  system.  However,  backups  occurred,  not  in  the  overall 
system  but  at  nodes  just  upstream  of  where  the  interruptions  took  place.  If  the 
density  was  large  enough  that  flow  reached  the  link  packet  transport  capacity 
without  sending  all  the  packets,  a  backup  occurred  and  grew  linearly  without 
bound  for  the  time  while  the  density  exceeded  what  could  be  sent.  After  the  sat¬ 
uration  conditions  are  removed  the  spike  decays.  Moreover,  for  interior  sources, 
the  combination  of  low  flux  and  upstream  input  could  also  lead  to  growth  at 
the  source  grid  point. 

The  work  (i)  derived  a  valid  model  for  a  discrete  buffer-link  model  of  com¬ 
munication  network  flow,  (ii)  demonstrated  how  the  FIFO  algorithm  is  imple¬ 
mented  in  such  a  model,  (iii)  derived  a  continuum  modes,  and  (iv)  provided  some 
solutions  to  the  model  demonstrating  how  it  can  describe  network  saturation 
and  interruptions. 
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